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Abstract 

We present algorithmic, complexity and implementation results for the problem of isolating the real 
roots of a univariate polynomial in Ba S L[y], where L = (Q(a) is a simple algebraic extension of the 
rational numbers. We revisit two approaches for the problem. In the first approach, using resultant 
computations, we perform a reduction to a polynomial with integer coefficients and we deduce a bound 
of C>b{N^°) for isolating the real roots of Ba, where N is an upper bound on all the quantities (degree 
and bitsize) of the input polynomials. In the second approach we isolate the real roots working directly 
on the polynomial of the input. We compute improved separation bounds for the roots and we prove that 
they are optimal, under mild assumptions. For isolating the real roots we consider a modified Sturm 
algorithm, and a modified version of dbscartbs' algorithm introduced by Sagraloff. For the former 
we prove a complexity bound of Ob{N^) and for the latter a bound of Ob{N^)- We implemented the 
algorithms in C as part of the core library of mathematica and we illustrate their efficiency over various 
data sets. Finally, we present complexity results for the general case of the first approach, where the 
coefficients belong to multiple extensions. 

Keywords real root isolation, algebraic polynomial, field extension, separation bounds, Sturm, Descartes' 
rule of sign 

1 Introduction 



Real root isolation is a very important problem in computational mathematics. Many algorithms are known 
for isolating the real roots of a polynomial with integer or rational coefficients that are either based solely 
on operations with rational numbers, js, 13, i^, 29 1 and references therein, or they follow a numerical, but 
certified approach, 2^, 33 1 and references therein. In this paper we consider a variation of the problem in 
which the coefficients of the polynomial are polynomial functions of a real algebraic number, that is they 
belong to a simple algebraic extension of the rationals. 



Problem 1. Let a be a real algebraic number with isolating interval representation a = {A, 3), where 
A = ET=o «^ J = [ai' 32], ai,2 e Q anddeg{A) = m andC{A) = t. Let = J27=o H") £ 
be square-free, wiiere hi{x) = Yl]'=o ^ ^[^li ^ {'^1,3) ^ '^i Vi < < i < d. What is the 

Boolean complexity of isolating the real roots of Be, ? 



Rump [sij, see also 30 1, presented an algorithm for the problem that is an extension of Collins and 
Loos jH] algorithm for integral polynomials. Johnson presented and compared various algorithms for 
Problem [TJ He considered a norm based algorithm that reduces the problem to root isolation of integral 
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polynomial (this is the approach that we consider in Sec. [3]) and extended three algorithms used for integral 
polynomials, i.e. Sturm (we present it in Sec. 14. 2p . the algorithm based on derivative sequence and Rolle's 
theorem and the algorithm based on Descartes' rule of sign [5] (we present a modified version in Sec. l4.3p . 
Johnson and Krandick [16] modified the latter and managed to replace exact arithmetic, when possible, 
with certified floating point operations; a novelty that speeds up considerably the computations. Along 
the same lines, Rouillier and Zimmermann [2^ presented an optimal in terms of memory used algorithm 
for integral polynomials that exploits adaptive multiprecision techniques that could be used for Problem [l] 
if we approximate the real algebraic number up to a sufficient precision. In a series of works 11, l3. 



a bitstream version of Descartes' algorithm was introduced. The coefficients of the input polynomial are 
considered to be real numbers that we can approximate up to arbitrary precision. We use the most recent 
version of this approach, which is due to Sagraloff [s^l, to tackle Problem [ij Last but not least, let us also 
mention the numerical algorithms due to Pan (i^l and Schonhage [s^ , that could be also used if approximate 
a in our problem up to a sufficient precision. 

Rioboo considered various symbolic algorithms for operations with real algebraic numbers, based on 
quasi Sylvester sequences. These algorithms could be used for Problem [U and they are closely connected 
with the Sturm algorithm that we present (Sec. 14.211 . However, we use different subalgorithms for sign 
evaluations and solving polynomials. The focus in ^] is on efficient implementation of the real closure in 

AXIOM. 

Problem [1] is closely related to real root isolation of of triangular systems and regular chains. In [i, 2C, 



38l , |39| algorithms and implementations are presented for isolating the real roots of triangular polynomial 



systems, based on interval arithmetic and the so-called sleeve polynomials. In the case of two variables the 
problem at study is similar to Problem [TJ In this line of research the coefficients of the algebraic polynomial 
are replaced with sufficiently refined intervals, hence obtaining upper and lower bounds (i.e. a sleeve) for the 
polynomial. Isolation is performed using evaluations and exclusion predicates that involve the non-vanishin: 
of the derivative. To our knowledge there is no complexity analysis of the algorithms. Neverthelss in 
evaluation bounds are presented, which are crucial for the termination of the algorithm, based on separation 
bounds of polynomial systems. However, the systems used for the bounds involve the derivative of the 
polynomial (this is needed for the exclusion criterion), which is not the case for our approach. In [3] the 
problem of real root isolation of 0-dim square-free regular chains is considered. A generalization of Vincent- 
Collins- Akritas (or Descartes) algorithm is used to isolate the real roots of of polynomials with real algebraic 
numbers as coefficients. This approach is similar to the direct strategy that we study. To our knowledge the 
authors do not present a complexity analysis since they focus on efficient algorithms and implementation in 

MAPLE. 

We revisit two approaches for isolating the real roots of a square-free polynomial with coefficients in a 
siraple algebraic extension of the rational numbers. The first, indirect, approach (Sec. [3]), already presented 



m 



17|, is to find a polynomial with integer coefficients which is zero at all roots of Be, isolate its real roots. 



and identify the intervals which contain the roots of Ba- We compute (aggregate) separation bounds for the 
resulting polynomial (Lem. [Tj), that are slightly better than the ones in [31], and prove that the complexity 
of the algorithm is Ob{N^'^), where N is an upper bound on all the quantities (degrees and bitsizes) of the 
input. The second approach (Sec. 14. ip is to isolate the roots of the input polynomial directly, using either 
the Sturm's algorithm or Sagraloff 's modified Descartes algorithm. We analyze the worst-case asymptotic 
complexity of the algorithms and we obtained a bound of Ob{N^) and Ob{N''), respectively. We obtain 
these complexity bounds by estimating improved separation bounds for the roots (Sec. 14.11 and Lem. [9]), 
that we also prove that they are optimal (Sec. 14. 4p . The bounds are better than the previously known ones 
l3, sB] by a factor of N . We empirically compare the performance of the indirect approach and the direct 



approach based on Sagraloff 's modified Descartes algorithm. The algorithms were implemented in C as part 
of the core library of mathematica, and we illustrate their behavior on various datasets (Sec. [5]). The 
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complexity bounds that we present are many factors better that the previously known ones. However, a 
fair and explicit comparison with the bounds in [17] is rather difficult, if possible at all, since, besides the 
improved separation bounds that we present, the complexity bounds of many sub-algorithms that are used 
have been dramatically improved over the last 20 years, and it is not clear how to take this into account in 
the comparison. 

Finally, we present a generalization of the first approach to the case where the input polynomials are 
univariate, but with coefficients that belong to multiple extensions (Sec. [6]). We derive (aggregate) separation 
bounds for this case (Lem. fT2|) and we sketch the overall complexity of the algorithm. The bounds are single 
exponential with respect to the number of extensions. 

Notation Ob means bit complexity and the Os-notation means that we are ignoring logarithmic factors. 
For A = X)f=i "^i^* £ deg(j4) denotes its degree. jC{A) denotes an upper bound on the bitsize of the 

coefficients of A, including a bit for the sign. For a G (Q, £ (a) > 1 is the maximum bitsize of the numerator 
and the denominator. 

If ad are the distinct, possible complex, roots of A, then Ai = |qi — a^l, where ac, is the roots 

closest to Qi. A{A) = mini Ai(j4) is the separation bound, that is the smallest distance between two (real 
or complex, depending on the context) roots of A. By S(j4) = — Y2^=i lg^i(^)i '^^ denote the numbers of 
bits needed to represent isolating rational numbers for all the roots of A. 

Given two polynomials, possible multivariate, / and g, then resx{f,g) denotes their resultant with 
respect to x. 



2 Preliminaries 

Real algebraic numbers are the real roots of univariate polynomials with integer coefficients; let their set 
be IRaig. We represent them in the so-called isolating interval representation. If a £ IRaig then the 
representation consists of a square-free polynomial with integer coefficients, A £ 2[a;], that has a as a real 
root, and an isolating interval with rational endpoints, J = [ai,a2], that contains a and no other root of the 
polynomial. We write a = {A, J). 

The following proposition provides various bounds for the roots of a univariate polynomial. Various 
versions of the proposition could be found in e.g. l3, 3^1 . We should mention that the constants that 
appear are not optimal. For multivariate bounds we refer to [l5|. 

Proposition 1. Let f be a univariate polynoniial of degree p. If ji are the distinct real roots of f, then it 
holds 

h\ < 2||/||oo < 2^+1 , (1) 
-lgA(/) < -^lg|3disc(/,,d)| + ^lg(p) + 

(p-l)lg||/,ed||2 (2) 

< 2plgp+pT , 
-^lgA,(/) < -^lg|disc(/,,rf)|-f ^ ~2^~^ + 

i 

{2p - l)lg\\fredh (3) 

< 3p^ + 3pT + Aplgp , 

where fred is the square-free part of f, and the second inequalities hold if we consider / £ Z[a;] and C (/) = r. 
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Proposition 2. Let f G 2[a;] have degree p and bitsize r. We compute the isolating interval representation 
of its real roots and their multiplicities in Ob{p^ + p'^t + P^t^) j^l. [il /. The endpoints of tJie isolating 
intervals have bitsize 0{p^ + Pt) and £-{jred) = 0{p + r), where fred is the square-free part of f. If 
N = max{f),r} then complexity bound for isolation becomes Ob{N^)- 

Proposition 3. 0, H] Given a real algebraic number a = (/, [a,b]), where C{a) = C{h) = 0{p^ -'tPt), 
and g £ 2[a;], such that de.g{g) = q, L {g) = a, we compute sigrL{g{a)) in bit complexity Ob{p<1 max{T, cr} + 
pmin{p,9}V). 

For the proofs of the following results the reader may refer to Let f,g& (Z[a;])[j/] such that deg^{f) = 
p, deg^{g) = q, degy{f),degy{g) < d, t = max(/: (/) , £ (5)). By SR(/,5; a) we denote the evaluation of 
the signed polynomial remainder sequence of / and g with respect to x over a, and by SRj(/, 5 ; a) the j-th 
element in this sequence. 

Proposition 4. We can compute res{f, g) w.r.t. x or y in OB{pqm.a.x{p,q}dT). 

Proposition 5. We compute SR(/, g ; a), where a £ Qu{oo} andC (a) = a, in OB{pq max{p, g}dmax{T, cr}). 
For the polynomials SRj(/,5; a) £ except for f,g, we have deg^(SRj (/, 5 ; a)) = 0{{p + q)d) and 

L{SIlj{f,g; a)) = 0{max{p,q}T + mm{p,q}a). 

3 Reduction to integer coefficients 
3.1 Some useful bounds 

The roots of Ba in Problem [l] are algebraic numbers, hence they are roots of a polynomial with integer 
coefficients. We estimate bounds on the degree and the bitsize of this polynomial, and we will use them to 
analyze the Boolean complexity of the real root isolation algorithm. 

Consider a real algebraic number a £ IRaig, in isolating interval representation a = {A, J), where A = 
g Oi 3 = [ai,a2], ai,2 £ Q and deg{A) = m and jC{A) = t. Since A is square-free, has m, possible 
complex, roots, say ai, a2, • • • , Q^m and after a (possible) reordering let a = Qi. 

Let Ba £ 2(a) [j/], be a univariate polynomial in y, with coefficients that are polynomials in a with integer 
coefficients. More formally, let Ba = Ya^o ^ii'^) ?/^ where bi{x) = Yl]Lo ^13 and rji < m, < i < d. The 
restriction r/i < m comes from the fact that 2(a) is a vector space of dimensior0 m and the elements of one 
of its bases are 1, a, ... , a'^^^. Finally, let C {Ba) = max^j £, {cij) = a. We assume that Ba is a square-free. 

Our goal is to isolate the real roots of Ba (Problem[l]). Since Ba has algebraic numbers as coefficients, 
its roots are algebraic numbers as well. Hence, there is a polynomial with integer coefficients that has as 
roots the roots of Ba, and possible other roots as well. To construct this polynomial, e.g. [Si J^^, E3i 
consider the following resultant w.r.t. x 

m 

R{y) = res,{B{x, y),A{x)) = {-ir^ Yi B{a,,y), (4) 

j=i 

where r] = max{r7i}, and B{x, y) £ y] is obtained from Ba after replacing all the occurrences of a with 
X. Interpreting the resultant using the Poisson formula, R{y) is the product of polynomials B{aj, y), where 
j ranges over all the roots of A. Our polynomial Ba £ 2(a)[j/] is the factor in this product for j = 1. Hence, 
R has all the roots that Ba has and maybe more. 

'if j4 is the minimal polynomial of a then the dimension is exactly m. In general it is not (computational) easy to compute 
the the minimal polynomial of a real algebraic number, thus we work with a square-free polynomial that has it as real root. 
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Remark 6. Notice that R{y) is not square-free in generd.!. For exampie consider the polynornid.! Ba — 
— a"^, where a is the positive root of A = — 3. In this case R{y) = TeSx{A{x) ^ B{x ^ y) = resa;(a;'^ — 
3,y^-x') = {y^-Sr. 

Using Prop. [m and by taking into account that rji < m, we get deg(i?) < mn and L{R) < m{T + a) + 
2mlg(4mn). We may also write deg(i?) = 0{mn) and C {R) = 0{m{a + r)). 

In order to construct an isolating interval representation for the real roots of Ba, we need a square-free 
polynomial. This polynomial, C{y) £ is a square factor of R{y), and so it holds deg(C) < mn and 

L (C) < m(r + a) + 3mlg(4mn), where the last inequality follows from Mignotte's bound [23j. 

Using the Prop, [fl we deduce the following lemma: 

Lemma 7. Let Ba be as in ProbJem[IJ The minimal polynomial, C £ 2[a;], of the, possible complex, roots 
of Ba, 7i, has degree < mn and bitsize < m{T + cr) + 3mlg(4mn)) or 0{m{T + a)). Moreover, it holds 

Y)i\ < 2'"*-'^"^'^-'^^'"'^*-^'""-' (5) 
-lgA(C) < m2n(r + a + 41g(4mn)) , (6) 
-^lgAj(C) < 3m^n{n + T + a + 6lg{Amn)) , (7) 



|7,| < 2°('"(^+'^» , (8) 
-lgA(C) = d{m^n{T + a)) , (9) 
S(C) = -^lgA,(C) = d{m^n{n + T + a)) . (10) 



3.2 The algorithm 



The indirect algorithm for Problem [2l follows closely the procedure described in the previous section to 
estimate the various bounds on the roots of Ba- First, we compute the univariate polynomial with integer 
coefficients, R, such that the set of its real roots includes those of Ba- We isolate the real roots of R and we 
identify which ones are roots of Ba- 

Let us present in details the three steps and their complexity. We compute R using resultant compu- 
tation, as presented in ([4]). For this we consider B as a bivariate polynomial in 2[a;,y] and we compute 
TeSx{B{x, y), A{x)), using Prop.H Since deg^(B) < m, degj^(S) = n, L {B) = a, deg^{A) = m, degy{A) = 
and C {A) = r, this computation costs O B{m^ n{a -\- r)), using Prop. ID 

Now we isolate the real roots of R- This can be done in OB{m*n^{mn^ -^-mnT+n^a+mr^ +ma^ +mTa)), 
by_Prop.[2l In the same complexity bound we can also compute the multiplicities of the real roots, if needed 

The rational numbers that isolate the real roots of R have bitsize bounded by 0{m^n{n -I- a -I- r)), which 
is also a bound on the bitsize of all of them, as Prop. [1] and Lem. [7] indicate. 

It is possible that R can have more roots that Ba, thus it remains to identify which real roots of R are 
roots of Ba- For sure all the real roots of Ba are roots of R- Consider a real root 7 of i? and its isolating 
interval [ci,C2]. If 7 is a root of Ba, then since Ba is square-free, by Rolle's theorem it must change signs if 
we evaluate it over the endpoints of the isolating interval of 7. Hence, in order to identify the real roots of 
R that are roots of Ba it suffices to compute the sign of Ba over all the endpoints of the isolating intervals. 

We can improve the step that avoids the non-relevant roots of R by applying the algorithm for chainging 



the ordering of a bivariate regular chain [2g]. However, currently, this step is not the bottleneck of the 
algorithm so we do not elaborate further. 
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Consider an isolating point of R, say cj G Q, of bitsize Sj. To compute the sign of the evaluation of Ba 
over it, we proceed as follows. First we perform the substitution y = Cj, and after clearing denominators, 
we get a number in 2[a], for which we want to compute its sign. This is equivalent to consider the 
univariate polynomial B{x,Cj) and to compute its sign if we evaluate it over the real algebraic number a. 
We have deg(S(a;, Cj)) = 0{m) and L {B{x, Cj)) = 0{a + nsj). Hence the sign evaluation costs Osirn^T + 
■m?a + m^nsj) using Prop. [31 Summing up over all Sj's, there are 0{mri), and taking into account that 

Sj = (D{m?n{a + r + n)) (Lem. [7]), we conclude that the overall complexity of identifying the real roots 
of Ba is OBim^n^ + m^nr -I- m^na + m^n?[a + r)). 

The overall complexity of the algorithm is dominated by that of real solving. We can state the following 
theorem: 

Theorem 8. The complexity of isolating the real roots of B ^ 2(a) [y] using the indirect method is 
Osi'm'^n^iiTin^ + mnr + li^a + mr^ + ma^ + mra)). If N = max{m, n, cr, r}, then the previous bounds 
become OsiN'^^)- 



If the polynomial Ba is not square-free then we can apply the algorithm of [37| to compute its square-free 
factorization and then we apply the previous algorithm either to the square-free part or to each polynomial of 
the square-free factorization. The complexity of the square-free factorization is O b{i^'^ n{a^ -\- t'^) -\- m'n? {a -\- 
r)), and does not dominate the aforementioned bound. 

4 Two direct approaches 

The computation of R, the polynomial with integer coefficients that has the real roots of B^ is a costly opera- 
tion that we usually want to avoid. If possible, we would like to try to solve the polynomial B^ directly, using 
one of the well-known subdivision algorithms, for example strum or DESCARTES and Bernstein, specially 
adopted to handle polynomials that have coefficients in an extension field. In practice, this is accomplished 
by obtaining, repeatedly improved, approximations of the real algebraic number a and subsequently apply 
DESCARTES or BERNSTEIN for polynomials with interval coefficients, e.g. 



The fact that we compute the roots using directly the representation of B^ allows us to avoid the com- 
plexity induced by the conjugates of a. This leads to improved separation bounds, and to faster algorithms 
for real root isolation. 

4.1 Separation bounds for Ba 

We compute various bounds on the roots of Ba based on the first inequalities of Prop.[lJ For this we need 
to compute a lower bound for |clisc(Sa)| and an upper bound for ||Sa||2. 

First we compute bounds on the coefficients on Ba- Let ai = a, 02,..., am be the roots of A. We 
consider the resultants 

f ^ \ 

V J 

It holds that 

m 

= al,Y[iz- H^k)) , 
fc=i 

where r] = max{r7i} < m. The roots of are the numbers foi(afc), where k runs over all the roots of A. We 
use Prop. [14] to bound the degree and bitsize of r^. The degree of ri is bounded by m and their coefiicient 
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are of bitsize < ma + mr + 5mlg(m). Using Cauchy's bound, we deduce 

2-m(T-mT-5mlg(m) ^ |fcj(a^)| < 2™'^+'"'^+^'" 's(™) 

for all i and fc. To bound |disc(BQ,)| we consider the identity 

disc(B„) =(-l)^"("-')^-^res^(B„,aB„(2/)/5j/) 

On (a) 

On (a) 

where the resultant, Rb G ^[q], can be computed as the determinant of the Sylvester matrix of Ba and 
dBa{y)/dy, evaluated over a. 

The Sylvester matrix is of size (2n— 1) x (2n— 1), the elements of which belong to 2[q]. The determinant 
consists of (2n — 1)! terms. Each term is a product of n — 1 polynomials in a of degree at most m — 1 and 
bitsize at most a, times a product of n polynomials in a of degree at most m — 1 and bitsize at most cr + lgn. 
The first product results a polynomial of degree (n — l)(m — 1) and bitsize (n — l)cr + (n — l)lgm. The 
second product results polynomials of degree n(m — 1) and bitsize nalgn + nlgm. Thus, any term in the 
determinant expansion is a polynomial in a of degree at most (2n — l)(m — 1), or 0(mn), and bitsize at 
most 4(2n — l)(jlg(mn) or 0{na). The determinant itself, is a polynomial in a of degree at most mn and 
of bitsize 4(2n - l)alg{mn) + (2n - 1) lg(2n - 1) < 5(2n - l)alg{mn) = d{na). 

To compute a bound on Rb{<^) we consider Rb as a polynomial in and we compute a bound on its 
evaluation over a. For this we use resultants. It holds 

m 

D = res,(^(i), y - Rb{x)) = a^^g(«-) JJiy - Rsia,)) . 

We notice that the roots of D £ Z[a;] are the evaluations of Rb over the roots of A. So it suffices to compute 
bounds on the roots of D. Using Prop. [14] we deduce that deg(D) < m and JC{D) < 13mncrlg(mn) + mriT 
or C (D) = 0{mn{a + r)). Using Cauchy bound, refer to Eq. ([l]), we conclude that 

Using this inequality and (fTl|l . we can bound |disc(BQ,)|, i.e. 

2" ISmncr Ig(mn) — 2mnT ^ |ciisc(.B )| ^ (^^Zrnna \g{rnri)-\-2rnnT (12) 

It remains to bound ||Ba||2. Using Eq. ()lip we get 

n 

\\Ba\\l < <(" + !) 22'"('^+^+5ig('")) . 

The previous discussion leads to the following lemma 
Lemma 9. Let Ba be as in Probleni[^ and be its roots. Then, it holds 

1^,1 < 2'"('"+'"+'^'em) _ (13) 
-lgA(B„) < 12mn(c7lg(mn) + T + 51gm) , (14) 
-^lgAi(S„) < 14mn(c7lg(mn) + T + 51gm) , (15) 

i 

or 

1^,1 < 2^('"(^+'"» , (16) 

-lgA(B„) = OimniT + a)) , (17) 

S(S<,) = -^lgA,(S„) = d{mn{T + a)) . (18) 
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4.2 The STURM algorithm 

Let us first study the STURM algorithm. We assume Be, as in Problem [l] to be square-free. To isolate the 
real roots of Ba using the STURM algorithm, we need to evaluate the Sturm sequence of B{a, y) and its 
derivative with respect to y, dB{a, y)/dy, over various rational numbers. For the various bounds needed we 
will use Lem. \9\ 

The number of steps that a subdivision-based algorithm, and hence STURM algorithm, performs to isolate 
the real roots of a polynomial depends on the separation bound. To be more specific, the number of steps, 
(#T), that STURM performs is (#T) < 2r -t- r Ig B + T,{Ba) [i, [13], where r is the number of real roots and 
B is an upper bound on the real roots. Using ()14p and (fT5|) we deduce that (#T) = 0{mn{T + a)). 

To complete the analysis of the algorithm it remains to compute the complexity of each step, i.e. the 
cost of evaluating the Sturm sequence over a rational number, of the worst possible bitsize. The latter is 
induced by the separation bound, and in our case is 0(mn(r -I- a)). 

We consider B as polynomial in 2[a;, y] and we evaluate the Sturm-Habicht sequence of B and over 
rational numbers of bitsize 0(mn(T -I- a)). The cost of this operation is OB{iTi^n'*'{T + a)) (Prop. [5]). 

It produces 0(n) polynomials in 2[a;], of degrees 0{mn) and bitsize 0{nT -\- na). For each polynomial 
we have to compute its sign if we evaluate it over a. Using Prop.[3]each sign evaluation costs OB{'m{m?' + 
n^)T -I- mri^a), and so the overall cost is O B{'mn{m? + n?)T + mn^a). If we multiply the latter bound with 
the number of steps, 0(mn(T -I- a)), we get the following theorem. 

Theorem 10. The complexity of isolating the real roots of B ^ 2(a) [j/] using the STURM algorithm is 

OBirri^n^im? + n'^){T^ + a^)), or dB{N^), where N = max{m, n, a, r}. 

4.3 A modified dbscartes algorithm 

We consider Sagraloff's modified version of Descartes' algorithm js^], that applies to polynomials with 
bitstream coefficients. We also refer the reader to 12, 2l|. 



As stated in ProblemfU let a be a real root of j4 = Yl^o '^i^^ ^ where 7^ and \ai\ < 2^ for 

< i < m, and let Ba = Yl7=o £ where = J2]'=o '^hj^^ ^ ^[^] i < and \c^j\ < 2"^ for 

< i < n and < j < rji, where we also assume that Ba is square-free. 

Let ^1, . . . be all (complex) roots of B, and Ai(Ba) := minj^il^j — By Theorem 19 of [33], the 
complexity of isolating real roots of Ba is 



dB{n{^{Ba) + nTB 



\2\ 



where 



bn{a) 



< 2^s and I](B„) = - Yli=i lg(Aj(S„)). From Lem. [9] we get that 

T,{Ba) < UTnn{T + alg{mn)) + nlgn = d{mn{T + a)) . (19) 



< 22ma+2mT+6mlg{m) f^j. ^jj j 



To compute a bound on tb, we use Eq. (fTTjl . It holds 
Hence, 

Tb < 2ma + 2mT + 6m lg(m) = d{m{a + r)) . (20) 
Finally, by combining (fT9|) and ([20|) , we deduce that the cost of isolating real roots of B is 

(5B(n(S(Sa) + nTB^) = dB{n{mnT + mnaf) 
= 6B{'m?n\a^ +T^)) . 

If N = max{m, n, a, r}, then the bound becomes C)b(#''). 
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It remains to estimate the cost of computing the successive approximations of 6i(a)/6„(a). The root 
isolation algorithm requires approximations of 6i(a)/6„(a) to accuracy of 0(S(Sa) + utb) bits after the 
binary point. Since |fei(a)/fe„(a)| < 2"^^, to approximate each fraction, for < i < n — 1, to accuracy L, it 
is sufficient to approximate bi{a), for < i < n, up to precision 0{L + tb). Hence, the algorithm requires 
approximation of bi{a), for < i < n, to precision 0(S(S) + utb)- By inequality (fTT|) . |fci(a)| > 2^"^^, and 
therefore it is sufficient to approximate hi{a) to accuracy + titb)- 

Approximation of Ci^a' to accuracy of L bits requires approximation of a to accuracy of L + Ig |cij | + 
lg(i) + [j - 1) Ig < L + a + lg(m) + (m - l)(r + 1) = d{L + a + mr) bits. Hence the accuracy of 
approximations of a required by the algorithm is 



Therefore, the bit complexity of computing the required approximations of bi{a) /bn{cx) is 

d{m^T^ + rii?mn{a + r)) = d{m^{mT^ + na + nr)) . 

Theorem 11. The bit complexity of isolating the real roots ofBa of Probleml^using the modified Descartes' 
algorithm in 0J is OB(m^n^((T^ + r^) + m^(mT^ + na + nr)), or Ob{N''), where N = max{m, n, a, r}. 

4.4 Almost tight separation bounds 

Let a be the root of A{x) = x"^ - ax"^^^ - 1, in (a, a + 1), for a > 3, m > 3. Then the Mignotte 
polynomial Ba{y) = j/" — 2(a'^j/ — 1)-^, where k = [(m — 1)/2J, has two roots in (1/a'^ — h, + h), where 

h = Q;-fc("+2)/2 < ^-(m-2)(n+2)/4^ 

If a < 2"^ and r = fi(lg(mn)), then — IgA(Ba) = fi(mnT), which matches the upper bound in (flSl) of 
Lem. [9l This quantity, fi(mnT), is also a tight lower bound for the number of steps that an subdivision 
based algorithm performs, following the arguments used in [l^l to prove a similar bound for polynomials 
with integer coefficients. 

5 Implementation and experiments 

We compare implementations of two methods of real root isolation for square-free polynomials over sim- 
ple algebraic extensions of rationals. The first method, ICF (for Integer Continued Fractions), performs 
reduction to integer coefficients described in Section 13.21 For isolating roots of polynomials with integer 
coefficients it uses the MATHBMATICA implementation of the Continued Fractions algorithm The second 
method, BMD (for Bitstream Modified Descartes), uses Sagraloff's modified version of Descartes' algorithm 
([i^l, see Section fOjl . The algorithm has been implemented in C as a part of the mathematica system. 

For the experiments we used a 64-bit Linux virtual machine with a 3 GHz Intel Core i7 processor and 6 
GB of RAM. The timings are in given seconds. Computations that did not finish in 10 hours of CPU time 
are reported as > 36000. 

Randomly generated polynomials For given values of m and n each instance was generated as follows. 
First, univariate polynomials of degree m with uniformly distributed random 10-bit integer coefficients were 
generated until an irreducible polynomial which had real roots was obtained. A real root r of the polynomial 
was randomly selected as the extension generator. Finally, a polynomial in Z[r, y] of degree n in y and degree 
m — 1 in r with 10-bit random integer coefficients was generated. The results of the experiments are given 
in Table [TJ Each timing is an average for 10 randomly generated problems. 



0{T,{Ba) + titb) = d{mn{a ^ t)) . 



By Lemmata 4.4, 4.5 and 4.11 of 




181] , the bit complexity of approximating a to accuracy L is 



0(mV2 + m?L) . 
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1 gontj tun 


— ^ 


— — 


— ^ 


— 


— 


10 


ICF 


0.003 


0.006 


0.013 


0.082 


0.820 


BMD 


0.002 


0.002 


0.003 


0.006 


0.019 


20 


ICF 


0.004 


0.010 


0.048 


1.49 


2.80 


BMD 


0.008 


0.008 


0.010 


0.017 


0.053 


50 


ICF 


0.014 


0.044 


0.271 


8.29 


20.5 


BMD 


0.046 


0.050 


0.061 


0.079 


0.213 


100 


ICF 


0.047 


0.173 


1.09 


33.1 


108 


BMD 


0.165 


0.206 


0.137 


0.246 


0.546 


200 


ICF 


0.144 


0.612 


4.90 


141 


626 


BMD 


0.746 


0.701 


1.00 


0.824 


1.55 


Table 1. Randomly generated polynomials 


n 


Algorithm 


m = 2 


m = 3 


m = b 


m = 10 


m = 20 


10 


ICF 


0.011 


0.008 


0.032 


0.208 


1.75 


BMD 


0.007 


0.007 


0.009 


0.010 


0.015 


20 


ICF 


0.019 


0.041 


0.193 


1.50 


13.9 


BMD 


0.075 


0.071 


0.080 


0.088 


0.106 


50 


ICF 


0.122 


0.270 


1.51 


25.8 


338 


BMD 


1.78 


1.63 


1.83 


1.90 


2.27 


100 


ICF 


0.834 


2.17 


16.1 


365 


10649 


BMD 


54.7 


51.3 


56.0 


74.7 


92.4 


200 


ICF 


7.53 


31.2 


246 


8186 


> 36000 


BMD 


2182 


3218 


3830 


4280 


4377 



Table 2. Generalized Laguerre polynomials 



Generalized Laguerre Polynomials This example compares the two root isolation methods for gener- 
alized Laguerre polynomials L"(a;), where a was chosen to be the smallest root of the Laguerre polynomial 
Lm{x). Note that L'^{x) has n positive roots for any positive a and Lm{x) has m positive roots, so this 
example maximizes the number of real roots of both the input polynomial with algebraic number coefficients 
and the polynomial with integer coefficients obtained by ICF. The results of the experiment are given in 
Table H 



Generalized Wilkinson Polynomials This example uses the following generalized Wilkinson polynomi- 
als := nfc=i(^ ~ where a is the smallest root of the Laguerre polynomial Lm{x). The timings 
are presented in Table [31 



Mignotte Polynomials The variant of Mignotte polynomials used in this example is given by Mn^a{x) ■= 
yn _ 2(a'=y - 1)2, where a is the root of Am{x) := - Zx^'^ - 1 in (3,4), m > 3 and A; = [(m - l)/2j 
(see Section 113}. The results of the experiment are given in Table HI 

The experiments suggest that for low degree extensions ICF is faster than BMD, but in all experiments 
as the degree of extension grows BMD becomes faster than ICF. Another fact worth noting is that ICF 
depends directly on the extension degree m, since it isolates roots of a polynomial of degree mn. On the 
other hand, the only part of BMD that depends directly on m is computing approximations of coefficients, 
which in practice seems to take a very small proportion of the running time. The main root isolation loop 
depends only on the geometry of roots, which depends on m only through the worst case lower bound on 
root separation. Indeed, in all examples the running time of ICF grows substantially with m, but the 
running time of BMD either grows at a much slower pace or, in case of generalized Wilkinson polynomials, 
it even decreases with m (because the smallest root a of Lm{x), and hence the root separation of W„,a(a;), 
increase with m) . The superiority of the direct approach was also observed in [17| . 
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1 g o n t) ti in 


— — 


— — 


— — 


— — 2^ 


_m 20_ 


10 


ICF 


0.017 


0.012 


0.035 


0.285 


2.09 


BMD 


0.015 


0.013 


0.011 


0.015 


0.008 


20 


ICF 


0.029 


0.069 


0.262 


2.23 


18.3 


BMD 


0.059 


0.052 


0.069 


0.039 


0.027 


50 


ICF 


0.137 


0.356 


2.04 


45.4 


429 


BMD 


1.84 


1.35 


1.29 


0.703 


0.561 


100 


ICF 


0.808 


2.84 


24.6 


674 


8039 


BMD 


47.0 


38.6 


32.0 


23.3 


8.38 


200 


ICF 


8.48 


35.1 


348 


11383 


> 36000 


BMD 


3605 


2566 


2176 


927 


565 



Table 3. Generalized Wilkinson polynomials 



n 


Algorithm 


m = 3 


m = 5 


m = 10 


m = 20 


10 


ICF 


0.003 


0.008 


0.049 


0.594 


BMD 


0.010 


0.006 


0.014 


0.036 


20 


ICF 


0.006 


0.027 


0.288 


8.83 


BMD 


0.015 


0.020 


0.049 


0.137 


50 


ICF 


0.041 


0.441 


12.2 


777 


BMD 


0.112 


0.147 


0.321 


0.854 


100 


ICF 


0.866 


11.6 


729 


28255 


BMD 


0.702 


0.868 


2.32 


5.99 


200 


ICF 


35.7 


684 


23503 


> 36000 


BMD 


3.12 


5.30 


13.8 


46.1 



Table 4. Mignotte polynomials 



6 Multiple extensions 

In this section we consider the problem of real root isolation of a polynomials with coefficients in multiple 
extensions. We tackle the problem using a reduction to a polynomial with integer coefficients. The technique 
could be considered as a generalization of the one presented in Sec. [H 

We use x** to denote the monomial x^^ ■■■x^, with e = (ei,...,e^) £ IN^. For a polynomial / = 
Yl^=i Cjyi^^ e {ei, . . . , Grn} C IN^ be the support of /; its Newton polytope Q is the convex hull of 

the support. By (#(5) we denote the integer points of the polytope Q, i.e. (#(5) = IQ fl 



Problem 2. Let aj, where 1 < j < i, be a real algebraic numbers. Their isolating interval representa- 
tion is aj = {Aj,3j), where Aj = YliLo°'i^)' = [aj,i,aj,2], ai,2 e <Q, deg(j4j) = m, and C{Aj) = r. 
Let 

n 

=^6.(ai,...,a^)2/' £ 

2=0 

be square-free, where foi(x) = Cij x"^ £ ^[x], jC (cij) < a, for < i < d, and for e = (ei, . . . , e^), it 
holds ej < 7] < m, What is the Boolean complexity of isolating the real roots of Ba ? 



We denote by ai the coefficients of A^, where 1 < i < i, and by c the coefficients of B. We compute 
separation bounds following the technique introduced [if 

We consider the zero dimensional polynomial system (5) : Ai{x.) = ■ ■ ■ = yi^(x) = j4^_|_i(x) = 0, where 
^fc(x) = J2iLo '^'=,1 = 0, 1 < < and Ai+i = B(x, y) = Yli=o xg) j/' = 0. We should mention 

that we make the assumption that B does not become identically zero when Qi, . . . , are replaced with 
some set of their conjugates (otherwise the resultant is zero). 

We hide variable j/, that is we consider (5) as an overdetermined system oi l+l equations in I variables. 
We consider the resultant, R, with respect to ii, . . . , i^, that is we eliminate these variables, and we obtain 
a polynomial R £ 2[ai, . . . , a^, c, y\. We interpret the resultant using the Poisson formula [7|, see also |23], 
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R{y) = res^{Ai,...,Ai,B) = Y[B{ai,i^, . . . ,at^i,,y) , 

and R{y) E • • • ,a£,c])[y]. Similar to the single extension case, Ba, is among the factors of R, hence 

it suffices to compute bounds for the roots of R{y). 

We consider i? as a univariate polynomial in y. The resultant is a homogeneous polynomial in the 
coefficients of (S), we refer to e.g. [3, 27 1 for more details and to [IB] for a similar application. To be more 



specific, the structure of the coefficients of R is 

where 1 < A; < M^+i = m^, and i is a number in {1, . . . , n}. The semantics of a^' are that it is a monomial 
in the coefficients of Ai of total degree M^. Similarly, c^f^+i-'^ stands for a monomial in the coefficients of B 
of total degree M^+i — k. Moreover, < lrim^~^ < i{m — l)m^~^ < im^. The degree of R with respect to 
y is at most n M^+i = nm^. 
Since |aij | < 2"^, it holds 

IgHla^r- <T^W . (21) 

Similarly, since \cij \ < 2"^, we get 

lg|c|M^+i-fc < ^(^i _ < _ (22) 



Finally, |gfc| < Ylitl (#Qi)'^' M\, where (#(5i) is the number of integer points of the Newton polytope 
of the polynomial Ai. We let Ai+i = B. It is (#Qj) = m + 1 for 1 < i < ^, so HLi i#Qi)^' < 
(m + iYirn-i)m'-^ < ^^m' _ (#Q^+i) < (£(m - 1) + n)^+i Hence, 

(#Q,)M^+> <{{i{m-l) + nY+' +1+1)""' 

and so for every k 

Iglefcl < igYli*Q^)''' < 2lm'\g[mnl) . (23) 

2=1 

By combining (|2ip. ([22]) and ([23|l we can bound the coefficients of R and its square-free factors. Using 
also Prop. [1] we get the following lemma. 

Lemma 12. Let Ba be as in Problein\^ The minimal polynomial, Ci of the, possible complex, roots of 
Ba, 7i, has degree < nm^ and bitsize < m^{Tp + cr -I- 3^1g(mn^)) or 0{tv/-{Pt + cr)). Moreover, it holds 

|7j| < 2"'(^^'^+'^+^^'s(""^)) , (24) 
-lgA(C^) < m^^nifT a ^ Al\g{mnl)) , (25) 
-^lgA,(C^) < 'm?^n{fT + G + n-\-Ql\g{mnl)) (26) 

i 

|7,| < 2^("^'(^'^+'^)) , (27) 

-lgA(Q) = dim^'nieT + a)) , (28) 

-^lgA,(C,) = d{m^^n{eT + a + n)) . (29) 
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Remark 13. To match exaclty the bounds derived in Lem. Q one should use for Mi the more accurate 
inequaUty Mi < l{m - l)m^"^ 

We can isolate the real roots of d in dsin^rn^^ + n'^m^^Ti^ + n*m*^a + +n^m^^T^l'^ + +n^m^^a^). That 
is we get a single exponential bound with respect to the number of the real algebraic numbers involved. 
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A A bound for the resultant 



Proposition 14. Let B = Y^i j Cijx'^y^ £ 2[a;,?/] of degree n with respect to y and of degree rj with 
respect to x, and of bitsize a. Let A = X)i^o '^i^^ ^ ^[^] degree m and bitsize r. The resultant 
of B and A with respect to x is univariate polynomial in y of degree at most mn and bitsize at most 
ma + rjT + m lg(n + 1) + (m + rj) lg(m + rj) or 0{ma + jyr). 

Proof: The proof follows closely the proof in Prop. 8.15] that provides a bound for general multivariate 
polynomials. We can compute the resultant of B{x,y) and A{x) with respect to x from the determinant of 
the Sylvester matrix, by considering them as univariate polynomial in x, with coefficients that are polynomial 
in y, which is 



bn 



Orj-1 



bo 



ao 



\ 



bo 



ao 



bo 



V 



x'^-'^B 



X 
X°B 

x^-^A 
x^-'^A 



^m— 1 



ao 



where b^ = J27=o Ci,fc2/'- 

The resultant is a factor of the determinant of the Sylvester matrix. The matrix is of size {rj+m) x {rj+m), 
hence the determinant consists of [r] + m)\ terms. Each term is a product of m univariate polynomials in 
y, of degree n and bitsize a, times the product of n numbers, of bitsize r. The first product results in 
polynomials in y of degree at most mn and bitsize at most ma + mlg(n + 1); since there are at most 
(n + I)'" terms with bitsize at most ma each. The second product results in numbers of bitsize at most 
r/T. Hence each term of the determinant is, in the worst case a univariate polynomial in y of degree m and 
bitsize ma + r/T + mlg(n + 1). We conclude that the resultant is of degree at most mn in y and of bitsize 
ma + rjT + m lg(n + 1) + (m + rj) lg(m + rj) or 0{ma + rjT). □ 
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